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I. INTRODUCTION 

The hydrodynamical description of relativistic particle 
collisions was first discussed more than 50 years ago by 
Landau jlj and nowadays it is frequently used in differ¬ 
ent versions for simulations of heavy ion collisions. Such 
a simulation basically includes three main stages. The 
initial stage, the fluid-dynamical stage and the so-called 
Freeze Out (FO) stage when the hydrodynamical descrip¬ 
tion breaks down. During this latter stage, the matter 
becomes dilute, cold, and non-interacting, the particles 
stream toward the detectors freely, their momentum dis¬ 
tribution freezes out. Thus, the freeze out stage is essen¬ 
tially the last part of a collision process and the source 
of the observables. 

The usual recipe is to assume the validity of hydro¬ 
dynamical treatment up to a sharp FO hypersurface, 
e.g. when the temperature reaches a certain value, Tpo- 
When we reach this hypersurface, all interactions cease 
and the distribution of particles can be calculated. 

In such a treatment FO is a discontinuity where the 
properties of the matter change suddenly across some 
hypersurface in space-time. The general theory of discon¬ 
tinuities in relativistic flow was first discussed by Taub 
|2|. That description can only be applied to discontinu¬ 
ities across propagating hypersurfaces, which have space¬ 
like, {da^da^ = — 1), normal vector. The discontinuities 
across hypersurfaces with time-like, (dcr^da^ = 1), nor¬ 
mal vector were considered unphvsical. The remedy for 
this came only 40 years later in Q, generalizing Taub’s 
approach for both time-like and space-like hypersurfaces. 
Consequently, it is possible to take into account conser¬ 
vation laws exactly across any surface of discontinuity in 
relativistic flow. 

As it was shown recently in Hi, the frequently used 
Cooper-Frye prescription )g| to calculate post FO parti¬ 
cle spectra gives correct results only for discontinuities 
across time-like normal vectors. The problem of negative 


contributions in the Cooper-Frye formula was healed by a 
simple cut-off, Q(p IJ da M ), proposed by Bugaev jij. How¬ 
ever, this formulation is still based on the existence of a 
sharp FO hypersurface, which is a strong idealization of 
a FO layer of finite thickness Jjj. Thus, by assuming an 
immediate sharp FO process, the questions of final state 
interactions and the departure from local equilibrium are 
left unjustified. 

The recent paper 0, Q formulates the freeze out prob¬ 
lem in the framework of kinetic transport theory. The dy¬ 
namical FO description has to be based on the Modified 
Boltzmann Transport Equation (MBTE), rather than 
on the commonly used Boltzmann Transport Equation 
(BTE). The MBTE abandons the local molecular chaos 
assumption and the requirement of smooth variation of 
the phase-space distribution, f(x,p ), in space-time. This 
modification of BTE, makes it even more difficult to solve 
the FO problem from first principles. Therefore, it is very 
important to build phenomenological models, which can 
explain the basic features of the FO process. 

The present paper aims to build such a simple phe¬ 
nomenological model. The kinetic approach presented, 
is applicable for FO in a layer of finite thickness with a 
space-like normal vector. It can be viewed as a continua¬ 
tion and generalization of Earn 0 The kinetic model 
for FO in time-like direction was discussed in a recent 
paper jT3j, however, the fully covariant model analysis 
and the treatment are presented in Cl¬ 
in present work we use stationary, one-dimensional FO 
models for the transparent presentation. Such models 
can be solved semianalytically, what allows us to trace 
the effects of different model components, assumptions 
and restrictions applied on the FO description. We do 
not aim to apply directly the results presented here to 
experimental heavy ion collision data, instead our pur¬ 
pose is to study qualitatively the basic features of the 
FO process. We want to demonstrate the applicability 
of the proposed covariant FO escape rate, and most im- 


2 


portantly, to see the consequences of finishing FO in a 
finite layer. Up to now, two extreme ways of describing 
FO were used: (i) FO on an infinitely narrow hypersur¬ 
face and (ii) infinitely long FO in volume emission type 
of model. To our knowledge this is the first attempt 
to, at least qualitatively, understand how FO in finite 
space-time domain can be simulated and what will be its 
outcome. In such stationary, one-dinrensional models the 
expansion cannot be realistically included, therefore it is 
ignored. 

In realistic simulations of high energy heavy-ion reac¬ 
tions the full 3D description of expanding and freezing 
out system should be included. This work is under ini¬ 
tial development. 

II. FREEZE OUT FROM KINETIC THEORY 

A kinetic theory describes the time evolution of a single 
particle distribution function, f(x,p) = f(t,x,p°, p), in 
the 6D phase-space. To describe freeze out in a kinetic 
model, we split the distribution function into two parts 

UHE0: 

f(x,p) = fi(x,p) +ff(x,p). (1) 

The free component, //, is the distribution of the frozen 
out particles, while /,;, is the distribution of the interact¬ 
ing particles. Initially, we have only the interacting part, 
then as a consequence of FO dynamics, /,; gradually dis¬ 
appears, while ff gradually builds up. In this paper we 
convert the description of the FO process from a sudden 
FO, i.e. on a sharp hypersurface, into a gradual FO, i.e. 
in some finite space-time domain. 

Freeze out is known to be a strongly directed pro¬ 
cess EzL where the particles are allowed to cross the 
FO layer only outwards, in the direction of the normal 
vector, da^, of the FO hypersurface. Many dynamical 
processes happen in a way, where the phenomenon prop¬ 
agates into some direction, such as detonations, defla¬ 
grations, shocks, condensation waves, etc. Basically, this 
means that the gradients of the described quantity (the 
distribution function in our case) in all perpendicular di¬ 
rections can be neglected compared to the gradient in the 
given direction da^, i.e. y/ « da^d^f. In such a situa¬ 
tion these can be effectively described as one-dimensional 
processes, and the space-time domain, where such a pro¬ 
cess takes place, can be viewed as a layer. 

Therefore, we develop a one-dimensional model for the 
FO process in a layer of finite thickness, L. We assume 
that the boundaries of this layer are approximately par¬ 
allel, and thus, the thickness of the layer does not vary 
much. This can be justified, for example, in the case 
when the system size is much larger than L. At the 
inside boundary of this layer there are only interacting 
particles, whereas at the outside boundary all particles 
are frozen out and no interacting particles remain. Note 
that the normal to the FO layer, da can be both space¬ 
like or time-like. 


The gradual FO model for the infinitely long one- 
dinrensional FO process was presented in recent works 
camoa We are going to build a similar model, but 
now we make sure that FO is completely finished within 
a finite layer. 

A. Freeze out in a finite layer 

In kinetic theory the interaction between particles is 
due to collisions. A quantitative characterization of col¬ 
lisions is given by the mean free path (m.f.p.), giving the 
average distance between collisions. The m.f.p., A m fp, is 
inversely proportional to the density, A m fp ~ 1 /n(x). If 
we have a finite FO layer, the interacting particles inside 
this domain must have a finite m.f.p. During the FO pro¬ 
cess, as the density of the interacting particles decreases, 
they are entering into a collisionless regime, where their 
final m.f.p., tends to infinity, or at least, gets much larger 
than the system size L. The realistic FO process for nu¬ 
cleons in a heavy ion collisions happens within a finite 
space-time FO domain, which has a thickness of a few 
initial mean free paths |l8l . Hence, one must realize that 
the FO process cannot be fully exploited by the means of 
the m.f.p. concept, since we have to describe a process 
where we have on average a few collisions per particle be¬ 
fore freeze out. Therefore, this type of processes should 
be analyzed by having also another characteristic length 
scale different from the m.f.p. In our case it should be 
related to the thickness, L, of the FO layer. 

Recent conjectures based on strong flow and relatively 
small dissipation find that the state where collective 
flow starts is strongly interacting and strongly correlated 
while the viscosity is not large This indicates a 

small m.f.p., in the interacting matter, while at the sur¬ 
face A m f p —> oo. Several indications point out that in 
high energy heavy ion reactions freeze out and hadroniza- 
tion happens simultaneously from a supercooled plasma 
[20l [2lL E3 . This could be modeled in a way that pre¬ 
hadron formation and clusterization starts gradually in 
the plasma, and this process is coupled to FO in a fi¬ 
nite layer. The FO is finished when the temperature of 
the interacting phase drops under a critical value and 
all quarks cluster into hadrons, which no longer collide. 
This is the possible qualitative scenario with well defined 
finite thickness L of the FO layer. 

Now, let us recall the equations describing the evolu¬ 
tion in the simple kinetic FO model [l(J, [lj, [lj] ■ Starting 
from a fully equilibrated Jiittner distribution, fj(p), i.e. 
fi(s = 0,p) = fj(p) and //(s = 0 ,p) = 0, the two com¬ 
ponents of the momentum distribution develop in the 
direction of the freeze out, i.e. along da according to 
the following differential equations: 

d s fi(s,p)ds = —fi(s,p)W esc (s,p ), 

d s ff{s,p)ds = +fi(s,p)W esc (s,p) , (2) 

where W esc (s,p) is the escape rate governing the FO de¬ 
velopment and s = x^da^. Here x^ is a 4-vector having 
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its origin 1 at the inner surface, Si, of the FO layer, see 
Fig. □ In order to obtain the probability to escape, for 
a particle passing from 0 till s, 7 ? esc (s), we have to inte¬ 
grate the escape rate along a trajectory crossing the FO 
layer: 


fi(s,p) = fj(p)exp(-J ds'W esc (s',p)j . (3) 


l-Ve. 


The definition for the escape probability was previously 
given in 0, in terms of collision or scattering rates, 
where the FO process was lasting infinitely long. In our 
finite layer FO description the quantity that defines the 
escape probability is the escape rate. 



FIG. 1: The picture of a gradual FO process within the finite 
FO layer, in x-direction, i.e. da M = (0,1,0,0). The particles 
are moving in different directions outwards, indicated by the 
angle 9 . The inside boundary of the FO layer, Si (thick line) 
indicates the points where the FO starts. This is the origin of 
the coordinate vector, Within the finite thickness of the 
FO layer, L, the density of the interacting particles gradually 
decreases (indicated by shading) and disappears at the outside 
boundary, S 2 (thin line) of the FO layer. 

To have a complete physical FO finished at a finite 
distance/time, we require: V esc —> 1, when s —> L. In 
usual cascade models the probability of collision never 
becomes exactly zero, and correspondingly V esc never be¬ 
comes exactly one, and the FO process lasts ad infinitum. 
This is due to the fact that the probability of collision 
is calculated based on the thermally averaged cross sec¬ 
tion, which does not vanish for thermal, e.g. Gaussian, 
momentum distributions. In reality the free or frozen 


out particles have no isotropic thermal distributions but 
these distributions can be anisotropic and strongly con¬ 
fined in the phase-space. This means that the collision 
probability can be exactly zero and FO may be completed 
in a finite space-time domain. 

It seems reasonable to parameterize the escape rate, 
which has dimension one over length, in terms of some 
characteristic FO length, A'(s,p), 

W esc {x ,l ,p ,l ,d(Ty) = * ^ 0(p s ) > ( 4 ) 

where the cut-off factor, O (p s ) = <d{jp^da^), forbids the 
FO of particles with momenta not pointing outward 0|. 
This FO parameter, A'(s,p), is not necessarily an aver¬ 
age distance in space or duration in time between two 
subsequent collisions, like the m.f.p. The m.f.p., tends to 
infinity as the density decreases, while the FO just be¬ 
comes faster in this limit. Actually, the FO scale behaves 
in the opposite way to the m.f.p. This can be seen, for ex¬ 
ample, in a simple, purely geometrical freeze out model, 
which takes into account the divergence of the flow in a 
3D expansion !23j. Both this and the phase transition 
or clusterization effect described at the beginning of this 
section, lead to a finite FO layer L , even if the m.f.p., 
A m f p = 1/ncj is still finite at the outer edge of this layer. 

We consider the thickness of the layer L to be the 
’’proper” thickness of the FO layer, because it depends 
only on invariant scalar matter properties like cross sec¬ 
tion, proper density, velocity divergence, phase transition 
or clusterization rates. These should be evaluated in the 
Local Rest frame (LR) of the matter, and since the layer 
is finite, around the middle of this layer. The proper 
thickness is analogous to the proper time, i.e. time mea¬ 
sured in the rest frame of the particle, hence the proper 
thickness is the thickness of the FO layer measured in the 
rest frame attached to the freeze out front, that is, the 
Rest Frame of the Front (RFF). Some of the parameters 
like the velocity divergence and the phase transition rate 
describe the dynamical changes in the layer, so these can 
determine the properties, e.g. the thickness, of the finite 
layer. However, calculating L from the above mentioned 
properties is beyond the scope of this paper, and L is 
treated as a parameter in the following. 

Let us consider the Rest Frame of the Front, where the 
normal vector of the front points either in time, t , or in 
space, x, direction, introducing the following notations 2 . 
Indeed, if da M is space-like the resulting equations can be 
transformed into a frame where the process is stationary 
(here da M = (0,1, 0,0) and correspondingly s = x), while 
in the case of a time-like normal vector the equations can 
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Time-like 

Space-like 

1 Any point of the inner surface, Si, can be considered as an origin, 

dcr M 

(1,0, 0,0) 

(0,1,0,0) 

since translations along Si do not change s, the projection of x ^ 

s = ( d(j p x p ) 

t 

X 

on the FO normal vector, da^ , as long as Si and S 2 are parallel, 

p s = ( dcrppP ) 

p° 

p x 

as assumed. Of course, this latter assumption can be justified 

ds = ( da p d p ) 

dt 

d x 

only locally, in some finite region, as it is clear from Fig. HI 

\'(s) 

r'(t) 

\'(x) 
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be transformed into a frame where the process is uniform 
and time-dependent ( da^ = (1,0,0, 0), s = t). For sake 
of transparency and simplicity we will perform calcula¬ 
tions only for FO in a finite layer with space-like normal 
vector in this paper, but many intermediate results can 
and will be obtained in Lorentz invariant way. 

Inside the FO layer particles are separated into still 
colliding or interacting and not-colliding or free particles. 
The probability not to collide with anything on the way 
out, should depend on the number of particles, which are 
in the way of a particle moving outward in the direction 
p/|p| across the FO layer of thickness L, see Fig. |T| If 
we follow a particle moving outward form the beginning, 
= 0), i.e. the inner surface of the FO layer, Si, to a 
position there is still a distance 

L — s 

COS dp 

ahead of us, where 6 p is the angle between the normal 
vector and p/|p|. As this remaining distance becomes 
smaller the probability to freeze out becomes larger, thus, 
we may assume that the escape rate is inversely propor¬ 
tional to some power, a, of this quantity 0. fill . 

Based on the above assumptions we write the escape 
rate as: 


W esc = j (c°s6> p ) a 0(p s ), (5) 

where this newly introduced parameter, A, is the initial, 
i.e. at Si, characteristic FO length of the interacting 
matter, A = X'(s = 0,cos$ p = 1). The power a is influ¬ 
encing the FO profile across the front. Indeed, calculat¬ 
ing the escape probability, V esc , eq. 01, with the escape 
rate, given by eq. 0 , we find 

\ x cosS p O(p s ) 

) 



in further calculations. It is easy to check that our escape 
rate, eq. ©, equals the earlier expression 


P„ c = ^e( P ') 

in L —> oo limit. Thus, the model discussed in this paper 
is a generalization of the models for infinitely long FO, 
described in mmoi, and allows us to study FO in a 
layer of finite thickness. 

The angular factor, cos0 p , maximized the FO proba¬ 
bility for those particles, which propagate in the direction 
closest to the normal of the FO layer. For the FO in time¬ 
like directions, studied in [j]|, the angular factor was 1. 
This factor, and correspondingly the escape rate, eq. 0 , 
are not covariant. Furthermore, this earlier formulation 
does not take into account either that the escape rate 
of particles should be proportional to the particle veloc¬ 
ity (the conventional non-relativistic limit of the collision 
rate contains the thermal average < av >). Let us con¬ 
sider the simplest situation, when the Rest Frame of the 
Front is the same as Rest Frame of the Gas (RFG), where 
the flow velocity is = (1,0,0, 0). If freeze out prop¬ 
agates in space-like direction, i.e. da M = (0,1,0, 0), as 
shown in Fig. ^ then cos6* p = p®/]p|. Therefore, a 
straightforward generalization of the escape rate, based 
on the above arguments, is 


cos0 p x 


v = 


IPI 

„0 


V ) 


( 6 ) 


where the r.h.s. of this equation is an invariant scalar in 
covariant form. Now, we assume that this simple gener¬ 
alization is valid for any space-like or time-like FO direc¬ 
tion, even when RFG and RFF are different Horn. 

Based on the above arguments, we can write the total 
escape rate from eq. © in a Lorentz invariant 3 form: 


W esc = 


A \L- 




e(p s ) , 


( 7 ) 


for a = 1, and 


Vesc = 1 - exp 




1- 


L — s 


for a / 1. Thus, we see that for different a-values we 
have different FO profiles: 

a = 1: power like FO, 

a > 1: fast, exponential like FO, 


a < 1: no complete FO within the finite layer, since 
V ea c does not tend to 1 as s approaches L. 

In papers liQ the authors were using a = 1, and 
were modeling FO in an infinite layer. In order to study 
the effects of FO within a finite space-time domain, we 
would like to compare the results of our calculations with 
those of earlier works, therefore we shall also take a = 1 


which now opens room for general study of FO in rela¬ 
tivistic flow in layers of any thickness. 

Former FO calculations in minis were always per¬ 
formed in RFF. Aiming for semianalytical results and 
transparent presentation, as well as in order to compare 
our results with former calculations, we will also study 
the system evolution in RFF, but now this is only our 
preference. In principle calculations can be performed in 
any reference frame. In more realistic many dimensional 
models, which will take into account the system expan¬ 
sion simultaneous with the gradual FO, it will be proba¬ 
bly more adequate to work in RFG or in Lab frame, and 
our invariant escape rate, eq. 0, can be directly used 
as a basic FO ingredient of such models. 


3 Now, it is important that L is defined as an invariant scalar, so 
Wesc is also an invariant scalar. 
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B. The Lorentz invariant escape rate 


A) dap = (1,0,0,0), leads to W(p) = 1, 


In this section let us study this new angular factor, in 
more detail. We will take the p-dependent part of the 
escape rate, eq. 0 , and denote it as: 

W(p) = P -^ e^da,) . (8) 

P^Up 

In RFG, where the flow velocity of the matter is = 
(1,0,0, 0)rfg by definition, W(p), is given as 

W{p) = P d 0 °~ M 0(p p dcr M ) , (9) 

P RFG 

and it is smoothly changing as the direction of the nor¬ 
mal vector changes in RFG. This will be discussed in 
more detail in the rest of this section. 

In the following, we will take different typical points 
of the FO hypersurface, A, B, C, D, E, F, see Fig. |3 
At these points, the normal vectors of the hypersurface, 
dap = (h, i,j, k)RFG, are given below. 



B) dap = 7<t(1, v„, 0,0), leads to W(p) = t'Ap +v«p ) , 

C) dap = 7 e (l + e, 1 - e,0,0), where = ( 4 e)“ 5 , 
e <C 1. This leads to W{p) = 7 £ 0 +p ~ P \ 

D) dafj = 7 c(1 - e, 1 + e, 0, 0), where = (4e) _ 3, 
eCl. This leads to W(p) = 7 e ^ p +p o ~ P ^ 
x 0 (7 e (p° + p x ) - 7 e e(p° - p x )), 

E) dap = 'y <7 (v cr , 1,0, 0), leads to 

w{p) = 7 - ( % p 0 °~ px) x ©( 7ff (v aP ° - P x )y 

F) da^ = (0,1,0,0), leads to W(p) = ^ x Q(p x ). 

The resulting phase-space escape rates are shown in Fig. 
£3 for the six cases described above. 

Similar calculations can be done in RFF, where da p = 
(1,0,0,0) for A, B, C and da M = (0,1, 0,0) for D, E, F, 
leading the following values for W(p): 

A) up = (1,0,0,0), leads to W(p) = 1, 

B) = 7<t(1, —Va, 0,0), leads to 

^ (p) — ’ 

c) up = 7e(l + e,-(l - e), 0,0), where 7 £ = ( 4 e)“ 5 , 
e < 1. This leads to W{p) = ± {pU _ pX) ^ t[pHpX) , 

D) Up = 7 e (1 + e ,-(1 - e), 0 , 0 ), where 7 e = ( 4 e)“ 5 , 
e < 1. This leads to W{p) = ± {pa _ pX) p ° t{p0+pX) 

x 0(p x ), 

E) Up = 7 cr(l, —v„, 0,0), leads to 
W (P) = y a (j/-v^) X W). 

F) Up = (1,0,0,0), leads tolf(p)=^x Q(p x ). 


FIG. 2: A simple FO hypersurface in the Rest Frame of the 
Gas (RFG: [t,x]), where = (1, 0, 0, 0)rfg, including time¬ 
like and space-like parts. The normal vector of the FO front, 
dap, is a time-like 4-vector at the time-like part and it is 
changing smoothly into a space-like 4-vector in the space-like 
part. On these two parts of the hypersurface, in the Rest 
Frame of the Front (RFF), dap points into the direction of 
the t' (z'j-axis respectively. At points A, B, C, D, E, F, we 
have different Rest Frame(s) of the Front, (RFF: [t’,x’]). 


To calculate the normal vector for different cases shown 
in Fig. |21 we simply make use of the Lorentz transforma¬ 
tion. The normal vector of the time-like part of the FO 
hypersurface may be defined as the local f'-axis, while 
the normal vector of the space-like part may be defined 
as the local aZ-axis. As dap is normalized to unity, its 
components may be interpreted in terms of 7 ^ and v a , 

where 7 ^ = . 1 = . So, we have: 

V 1 ~ V S 


For these cases, A, B, C, D, E, F, in RFF the resulting 
phase-space escape rates are shown in Fig. Q] 

Figs. [3]andQ]show that the momentum dependence of 
the escape rate, uniform in point A, becomes different at 
different points of the FO hypersurface, but this change 
is continuous, when we are crossing the light cone, from 
point C to point D. Although in RFF, Fig. QJ it seems 
that there is a principal difference between space-like and 
time-like FO directions, due to the cut-off Q(p^dap) func¬ 
tion, but this is only the consequence of the chosen refer¬ 
ence frame, i.e. RFF is defined in a way to stress the dif¬ 
ference between these two cases, since going from C to D, 
the normal vector has a jump, i.e. dap = ( 1 , 0 , 0 , 0 ) goes 
over to da p = (0, 1 , 0, 0). Nevertheless, W(jp) is a contin¬ 
uous function as we change dap, and in other frames, for 
example in RFG, Fig. [31 we can see this clearly. 
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Px 


D E F 





Px 


FIG. 3: The contour plots of the momentum dependent part 
of the escape rate, W(jp ), as in eq. itift . presented in six sub¬ 
plots at different points of the FO hypersurface. All plots 
are in RFG. For region A: da M = (1,0,0,0), and W(p) 
is one uniformly, for B: dcr M = 7 CT (1, 0.5, 0, 0), C: da M = 
7e(1.01, 0.99, 0, 0), (first row), D: da M = 7 e (0.99, 1.01, 0, 0), 
E: da^ = 7 CT (0.5, 1, 0, 0), F: da^ = (0,1, 0,0), (second row). 
The momenta are in units of particle mass, [m]. 
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cT 0 
-2 


-2 0 2 




Px 
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FIG. 4: The contour plots of the momentum dependent 
part of the escape rate, W(p), in RFF. For region A: 
Ufj, = (1,0,0,0), and W(p) is one uniformly, for B: u M = 
7 ct (1, — 0.5, 0,0), C: = 7 e (1.01, —0.99, 0, 0), (first row), 

D\ un = 7 e (1.01, — 0.99, 0, 0), E: = 7 CT (1, —0.5, 0, 0), F: 

= (1,0, 0,0), (second row). In cases D, E, F the escape 
rate vanishes for momenta with p x < 0. The momenta are in 
units of particle mass, [m]. 


C. The updated simple kinetic model 


Now, using the new invariant escape rate, ecr J7J, we 
can generalize the simple model presented in [iJJ, UH E3 ■ 


i.e. eqs. ©, for a finite space-time FO layer: 


dsfids = - 

d s ffds = 


L — sj \p^u 


L 


B(p s )/, 


ds 

T’ 


L — s) \p^u 


0(p s )/, ^ , (10) 


Solving the first equation we find for the interacting com¬ 
ponent: 


fi{s,p) = fj{p) 

s—>L 


L — s 


L(-t 

\\pP 


e(p s ) 


fj(p)@(~P s )• 


( 11 ) 


Now, inserting this result into the second differential 
equation, from eqs. 113, we obtain the FO solution, 
which describes the momentum distribution of the frozen 
out particles: 


ff(s,p) = hip ) 


1 - 




/j(t)©(t s ) ■ 


( 12 ) 


As s tends to L , i.e. to the outer boundary of the FO 
layer, this distribution, depending on the direction of 
the normal vector (space-like) or time-like will tend to 
the (cut) Juttner distribution, fj(p). This means that 
(part of) the original Juttner distribution survives even 
when we reach the outer boundary of the FO surface. To 
remedy this highly unrealistic result, in Em El El, 
rethermalization in the interacting component was taken 
into account via the relaxation time approximation, i.e. 
we insert into the equation for the interacting compo¬ 
nent a new term, which describes that the interacting 
component approaches some equilibrated (Jiittner) dis¬ 
tribution, f eq (s ), with a relaxation length, Aq: 


dsfids = 


+ 


dsffds = 


P 


L — s 
[feq{s) - h] 


p^U/j. 

ds 

Aq 


e(p*)/, * 




e(p’)/, 


ds 

T' 


(13) 


(14) 


Let us concentrate on the equation for the interacting 
component. Here the first term from eq. 113. related to 
FO, moves the distribution out of the equilibrium, and 
decreases the energy-momentum density and baryon den¬ 
sity of the interacting particles. The second term from eq. 
113, changes the distribution in the direction of the ther- 
malization, while it does not effect the conserved quan¬ 
tities. The relative strength of the FO and rethermal¬ 
ization processes is determined by the two characteristic 
lengths, A and Ao- 

In general the evolution of the interacting component 
can be solved numerically or semianalytically, at every 
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step of the integration. Then, the change of conserved 
quantities due to FO should be evaluated using the actual 
distribution, /,;(s,p) at the corresponding point s. For 
the purpose of this work, namely for the qualitative study 
of the FO features, it is enough to use an approximate so¬ 
lution, similarly as it was done in [ill LL2 j ll5■ This would 
also allow us to make a direct comparison with results 
of these older calculations. Thus, the evaluation of the 
change of the conserved quantities is done analytically, 
i.e. /i(s,p) is approximated with an equilibrium distribu¬ 
tion function f e q{s) with parameters, T(s), n(s ), u M (s). 

This approximation is based on the fact that in most 
physical situations the overall number of particle colli¬ 
sions vastly exceeds the number of those collisions, after 
which a particle leaves the system or freezes out. This 
allows us to take that rethermalization 4 happens faster 
than the freeze out, i.e. that Ao < A or Ao <C A. Of 
course, this argument is true only at the beginning of 
the FO process, when the density of the interacting par¬ 
ticles is still large. When s is close to L, i.e. near the 
outer hypersurface, the first term in eq. m becomes 
more important than the rethermalization term because 
of its denominator, but as we shall see in the results sec¬ 
tion, particles freeze out exponentially fast and for large 
s, when say 99% of the matter is frozen out, the error we 
introduce with our approximate solution can not really 
affect the physical situation. 



FIG. 5: Numerical solutions of eq. m for different Ao’s. 

_ 2 

This solution was obtained using, f eq = e 3 , test function, 
for L = 10A. The results show that fi(s) approaches f eq (s) 
when A 0 A. 

For illustration, let us take a test function, / eq (s) = 
e~ s , (we ignore p dependence for the moment), which 
is a smoothly and fastly decreasing function 5 . In Fig. [3 


4 The words ’’immediate rethermalization” used in a few earlier 
publications, were badly chosen, misleading and inappropriate. 

5 In real calculations the s-dependence of f eq (s,p) is calculated 
from the energy, momentum and baryon charge loss of the in- 


we show the numerical solutions for the interacting com¬ 
ponent for Ao = A and Ao = 0.1A. The results show 
that for the latter case we can safely take for the fi the 

approximate solution 0111113 : 

fi(s ) = feq(s) ■ (15) 

D. Conservation Laws 

The goal of the freeze out calculations is to find the 
final post FO momentum distribution, and then the cor¬ 
responding quantities defined through it, starting from 
the initial pre FO distribution. On the pre FO side we 
can have equilibrated matter or gas. Its local rest frame 
defines the RFG, see Fig. El We can also define the ref¬ 
erence frame, which is attached to the freeze out front, 
namely the RFF, see Fig. [3 These choices are usually 
advantageous, but other choices are also possible. 



FIG. 6: The orientation of the freeze out front in RFG is 
given as, da M = 7 a ( ua , 1,0, 0)rfg- 

Furthermore, the conservation laws and the nonde¬ 
creasing entropy condition must be satisfied [lj} : 

[mda fi ) = 0, [T^da^]=0, [S^da^ > 0, (16) 

where [A] = A—Aq. The pre FO side baryon and entropy 
currents and the energy-momentum tensor are denoted 
by N[f, Sq,Tq' / , while the post FO quantities are denoted 
by N^, 

The change of conserved quantities caused by the parti¬ 
cle transfer from the interacting matter to the free matter 
can be obtained in the following way. For the conserved 
particle 4-current we have 

dN^ = dN? + dNf = 0 =7 dN? = -dNf. (17) 


teracting component, where these losses are determined by the 
momentum dependence of the escape rate and the actual shape 
of fi(s,p ), as discussed here. 

















Then, using the kinetic definition of the particle current, 
together with eqs. mm, we obtain 


dN?(s ) 



Similarly the change in the energy-momentum is 


dTr(s) 


ds 

X 



' V s 
.Pp u i 


© (p s ) 


feq(s) ■ 


The parameters of the equilibrium (Jiittner) distribution, 
feq(s), have to be recalculated after each step, ds , from 
the conservation laws as in fU HHH- 
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and for the change of energy density we have 


dei{s) = u^i(s) dT^ v (s) u v ,i(s). (21) 


The change of the temperature of interacting component 
can be found from this last equation, eq. m, and from 
the Equation of State (EoS). This closes our system of 
equations. 

If we fix the FO direction to the ^-direction, then eqs. 
(ITTO can be rewritten as: 


dKM = -XtXi / $v 


and 


A \L — x 


p cos 0 V 


pj 


( 22 ) 


dTr(x) = - 


i{p° — j u p cos d p ) 
r d 3 p 


0 (cos 8 P ) 


feq{x,p) . 


dx ( L 


A \L — xj J p° 
p cos 0 P 

7 (p° — jup cos 0 P ) 


P> v 


(23) 


0(cos Op) 


feq{x,p) 


where, the four momentum of particles is p ^ = (p°,p), 
p= |p|, p x = p cos Op, the flow velocity of the interacting 
matter is uf = 7 ( 1 , v, 0 , 0 ), 7 = -y==?, u = |u| and 
j = sign{v) . 


FIG. 7: The orientation of the freeze out front in RFF, i.e. 

in its own rest frame, is do> = (0,1, 0, 0 )rff- In this frame Changes of the conserved current and 

the gas has nonvanishing velocity in general. energy-momentum tensor 


The change of flow velocity, duP{s), can be calculated 
using Eckart’s or Landau’s definition of the flow, i.e. from 
dNd(s) or dTd v {s) correspondingly. Then the change of 
conserved particle density is given by 

drii (s) = u?(s) dN z a ( s ), (20) 


In this section we show new analytical results for the 
changes of the conserved particle current and energy- 
momentum tensor. The formulae are analogous to those 
from ref. dam, but now they are calculated with the 
Lorentz invariant angular factor from eq. 0 . We show 
results for both massive and massless particles. 


dN?(x) = 


dx 


A \L — x J 4m 2 7 2 


7Gq (m) + 


& 3 T( 0 , b) 


— 27 bu 


- Uu 


(1 +j)K 0 (a) - /C 0 (a, b) 


+ 


37 2 

it 2 7 2 (l + b)(u 2 — 3) — 


(1 +j)K 1 (a) - /Ci (a, b) 

Mb)]S 


37 2 


dx / 


A V L — x J 4 


(3 — i>)(l + v) 3 _ 2 


m —0 
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dN?(x) = 


dN?{x) dx ( L 


- 2j^bu 2 (1 +j)K 1 (a) - JCi(a,b) 


ju A \L — x) 4 w 2 7 2 

\ 

jb 2 u 2 \{l + j)K 0 (a,b) - JC 0 (a,b)] - ju 7 2 (l + u 2 )(l + b)e~ b 


m—0 


dx f L An 


A \L — x ) 4 


2(1 + v) 3 2 


7 


<«*(*) = 


dx ( L A nT j b 4 T{0, b) 


A \L — x) 4m 2 7 2 | 47 

— it6 3 [(1 + j)A'i(a) —/Ci(a,6)j +7 3 | — A(b)(3u 2 + 1) + 


3 - lG 2 (to) - u'yb (u + 3) |^(1 + j)K 2 (a) - IC 2 (a, &) 

B(b ) b 2 (3 + b) 


27 ® 


37 ® 


i,V(1 . + l,) + |i,V(„ 2 - 3) + a 2 

o 


m=0 




dx ( L \ nT 




A \L — x J 4 


(1+x) 4 


(6 — 4z> + t; 2 ) 7 3 


(24) 


_ dTfV) dx [ L \ nT , 2 2 


<HT(*) = 


ju X \L — x J 4 w 2 7 2 

V 

- j&(6V - 2) f(l + j)A 1 (a) - AC!(a, 6) 


-Jib (3 u + 1)\(1 + j)K 2 (a) - JC 2 (a,b) 


JU'y 


+ jab 
4 b 3 


(l+j)K 0 (a) - AC 0 (a,6) 


-y 2 A(b)(u 2 + 3) + (1 + 6 )(u 2 - 3) + —u 2 j 2 - b 2 

O 


-6 


ra =0 


dx f L \ nT 


A \L — x J 4 


(1 + f ) 4 


(4 - x) 7 3 


dT 0x (r) T 
dT xx (x) = ( ’ - 2 — 

ju 7 J u 


dNf - 


dNf (x) 


ju 


dx ( L 


nT 


A \L — x J 4 w 2 7 2 


- 7 it 6 2 (3 + w 2 ) (1 + j)K 2 (a) - £ 2 ( 0 , 6 ) 


A{b){3u 2 + 1) + -b 3 u 4 (u 2 — 3) + a 


-b 


m =0 


dx ( L \ nT 


A \L — x J 4 


3(1 + *Q 4 3 

o ^ 


dT[ v (x) = - 


dT xx {x) dx ( L 


nT 


A \L — x) Su 2 y 


- 7«& 2 (3 + u 2 ) (1 + j)K 2 (a) - IC 2 (a,b) 


A(b){3u 2 + 1) — 


35(6) {b 2 + 1)(1 + b) 2 


27 ® 


+ 


27 ® 


+ ^-b 3 u 4 (u 2 - 3) + a 2 

O 


- b b 4 no,b) 


e — 


47 s 


+ lH 2 ( to ) 


ra =0 


dx 


?zT 


A \L — x J 


7(1 +^) 3 
4v 2 


(-4 + 12v - 9v 2 + 3v 3 ) 


and 


dTr(x) = dTf v (x ), 


( 25 ) 
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where, a = b = try, A(b) = (2 + 2 b + b 2 )e b , B(b) = 
|(6 + 66 + 36 2 + b 3 )e~ b and n = 4TrT 3 a 2 K 2 {a) g 
is the particle density, g is the degeneracy factor while, 
G~{m), H~(m ), /C(a, 6), K(a), T(0, 6) are defined in Ap¬ 
pendix A. Note that the ^-dependent factor L/{L — x ) 
is just a multiplier in these calculations, and tends to 
unity if we are dealing with an infinitely long FO, as in 

E3EM3- 

III. RESULTS AND DISCUSSION 

In this section we calculate the post FO distributions 
and compare the results to former calculations presented 
in |1(J, |nl [12}. The effect of two main differences due to 
the new Lorentz invariant escape rate, eq. 0 , is to be 
checked: 

- the infinite, (00), FO layer or finite, (L), FO layer, 


Pl- Next, we are using the simple angular factor, but 
in this case inside a finite FO layer, L = 10A, 

Pl = { t - b) £ 5 r Ee( ' ,cose -)' 

The results are shown in Figs. m mm 

Woo'. Then, we are dealing with the new Lorentz invari¬ 
ant angular factor in the escape rate and with an 
infinite FO length, x max = 300A, 

Woo = y f pCb>s0p \ Q(pcos 8 p ). 

A V P ^ u m / 

The results are shown in Figs. 0031130. 

11+: Finally we present the primary results of this pa¬ 
per, using both our new improvements, i.e. the 
covariant escape rate of eq. 0 , 


- the simple angular escape rate, P , or the covariant 
escape rate, W. 


W L = 


A V L — x 


p cos 6 P 

P^Uft 


8(pcos 9 P ). 


We performed calculations for a baryonfree massless gas, 

where we have used a simple EoS, e = (JsbT a , where 
2 

&SB = 70- The change of temperature is calculated 
based on this EoS and eq. <ED- There are no conserved 
charges in our system, consequently we use Landau’s def¬ 
inition of the flow velocity [ITj : 


du i,Landan( X ) 


A^{x) dTj^ a (x) u°(x) 
ei(x) + Pi(x) 


(26) 


where A^ v (x) = g— u^(x) u^(x) is a projector to the 
plane orthogonal to u^(x), while ei(x ) and Pi(x ) are 
the local energy density and pressure of the interacting 
component, i.e. T? v (x) = ( ei(x ) + Pi(x))u^(x)v%{x) — 
Pi(x)g^. A detailed treatment of Eckart’s flow velocity 
can be found in mm. 

For such a system we finally obtain the following set 
of differential equations: 


dlnT = 
dv = 


7 


4<t S bT 4 

3 

4 usbT 4 


dT°° - 2 vdT? x + v 2 dT, xx 


- vdTr + (1 + v 2 )dT" x - vdT) 


(27) 


We will present the results for four different cases: 

Poo- We use the simple, but relativistically not invariant 
angular factor, cos0 p , in the escape rate, 


Poo = C °^ p 0(pcos dp) . 

The system is characterized by an infinite FO 
length (up to x ma x = 300A in calculations). The 
results are shown in Figs. 0 mm- This is the 
same model as in [T(J 10, ll3L 


The results are shown in Figs. 0 ED HMD- 

In the case of FO in the infinite layer the factor, L/(L—x), 
was replaced by 1. We presented the situation at a dis¬ 
tance of Xm ax = 300A, where the amount of still inter¬ 
acting particles is negligible. 

For the particular cases when we are dealing with an 
infinite FO, i.e. Poo and Woo, or with finite layer FO, 
i.e. Pl and Wl, the results are plotted together. Thus, 
in one figure the focus is on the consequences caused by 
the different angular factors. 

Thin lines always denote the cases with simple rela¬ 
tivistically not invariant angular factor, these correspond 
to Poo and Pl . Thick lines always correspond to cases 
with covariant angular factor, Woo and Wl- 
All the figures are presented in the RFF. 


A. The evolution of temperature of the interacting 
component 

The first set of figures, Figs. 00 shows the evolution 
of temperature of the interacting component, in fact the 
gradual cooling of the interacting matter, for the differ¬ 
ent cases, Poo, Pl, W x and Wl- 

First, on all figures matter with larger (positive) flow 
velocity, vq, cools faster. This is caused by the momen¬ 
tum dependence of the escape rate, which basically tells 
that faster particle in the FO direction, will freeze out 
faster. Thus, the remaining interacting component cools 
down, since the most energetic particles freeze out more 
often than the slow ones. Of course, for larger initial flow 
velocity, vq, in the FO direction, there are more particles 
moving in the FO direction with higher momenta in av¬ 
erage, than for a smaller flow velocity. 

Now, comparing Fig. 0with Fig. 0 we can see the 
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difference between finite and infinite FO dynamics. In a 
finite layer the cooling of interacting matter goes increas¬ 
ingly faster as FO proceeds, while for FO in infinite layer 
the cooling gradually slows down as x increases. The rea¬ 
son is the factor, L/(L — x), which speeds up FO as L — x 
decreases, and forces it to be completed within L. 



FIG. 8: The temperature of the interacting component in 
RFF for a baryonfree massless gas, calculated with the two 
escape rates: P a0 (thin lines) and Woo (thick lines), for an 
infinitely long FO, (xmax = 300A). The initial temperature is 
To = 170 MeV, Vo is the initial velocity in RFF. 



FIG. 9: The temperature of the interacting component in 
RFF, calculated with the two escape rates: Pl (thin lines) 
and Wl (thick lines), for a finite, (L = 10A), FO layer. The 
initial temperature is To = 170 MeV, Vo is the initial velocity 
in RFF. 

The difference between P and W escape rates comes 
from the denominator, which is p° in case of P and p^u^ 
in case of W. This difference leads to a stronger cooling 
for the escape rate P which is bigger than W, if vq ^ 0. 
This can be seen well at later stages of infinitely long FO, 
Fig. |dl particularly for the positive initial flow velocity. 
In all other cases the difference between old and new an¬ 


gular factors is insignificant, what supports our ’’naive” 
generalization of the angular factor. 


B. The evolution of common flow velocity of the 
interacting component 

The second set of figures, Figs. E3 El shows the evo¬ 
lution of the flow velocity of the interacting component. 

In both cases the flow velocity of the interacting com¬ 
ponent tends to —1, because the FO points to the positive 
direction and particles with positive momenta freeze out. 
Thus, the mean momentum of the rest must become neg¬ 
ative. 



FIG. 10: The flow velocity of the interacting component in 
RFF for a baryonfree massless gas, calculated with the two 
escape rates: Poo (thin lines) and Woo (thick lines), for an 
infinitely long FO, (i mas = 300A). The initial temperature is 
To = 170 MeV, vo is the initial velocity in RFF. 



FIG. 11: The flow velocity of the interacting component in 
RFF, calculated with the two escape rates: Pl (thin lines) 
and Wl (thick lines), for a finite, (L = 10A), FO layer. The 
initial temperature is To = 170 MeV, no is the initial velocity 
in RFF. 
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Comparing Fig. m with Fig. m we can see again 
that in a finite layer the flow velocity decreases faster 
and faster as FO proceeds, while for FO in infinite layer 
the velocity change gradually slows down as x increases. 
The reason is the L/(L — x) factor, as discussed above. 

The difference between the evolution of the flow ve¬ 
locity, due to the different angular factors, is again not 
significant, supporting its generalization. 

C. The evolution of the transverse momentum and 
contour plots of the post FO distribution 

The next set of figures, Figs. El El shows the evo¬ 
lution of the transverse momentum distribution, while 
Figs. El El present the contour plots of the post FO 
momentum distribution, for W^, and Wl- We have pre¬ 
sented a one-dimensional model here, but we assume that 
it is applicable for the direction transverse to the beam 
in heavy ion experiments. The presented plots should 
be qualitatively compared to the transverse momentum 
distributions of measured pions. 



FIG. 12: The local transverse momentum (here p x ) distribu¬ 
tion for a baryonfree massless gas at ( p y = 0), calculated with 
the two escape rates: Poo (thin lines) and Woo (thick lines), 
for an infinitely long FO, ( x ma x = 300A). The initial parame¬ 
ters are no = 0 and To = 170 MeV. The transverse momentum 
spectrum is obviously curved due to the freeze out process. 
The slope of the transverse momentum distribution increases 
as we are approaching infinity. 

What we see is that all the final post FO momentum 
distributions are essentially the same. This is very im¬ 
portant outcome from our analysis, which we will discuss 
below. Also, one can see that resulting post FO distribu¬ 
tions are non-thermal distributions, as it has been shown 
already in mmm, they strongly deviate from expo¬ 
nential form in the low momentum region. The increase 
in the final FO spectra over the thermal distribution for 
low momenta is connected to the fact that at late stages 
of the FO process, the interacting component is cold and 
its flow velocity is negative. So, it contributes only to the 


low momentum region of the post FO spectra. 


IX 10 X 300 X 



FIG. 13: The post FO distribution, ff(x, p), in RFF. The 
calculations were done with the Lorentz invariant escape rate, 
Woo, for an infinitely long FO, (Xmax = 300A). The subplots 
correspond to x = 1A, 10A, 300A respectively. The initial 
parameters are vq = 0 and To = 170 MeV. Contour lines are 
given at values represented on the figure. The maximum is 
increasing with x as indicated in Fig. El The distribution 
is asymmetric and elongated in the FO direction. This may 
lead to a large-pt enhancement, compared to the usual Juttner 
assumption used in many earlier calculations as a post freeze 
out distribution. Note that ff(x, p) does not tend to the cut 
Juttner distribution even at very large x. 



FIG. 14: The local transverse momentum (here p x ) distribu¬ 
tion for a baryonfree massless gas at ( p y = 0), calculated with 
the two escape rates: Pl (thin lines) and Wl (thick lines) for 
a finite, (L = 10A), FO layer. The initial parameters are 
no = 0 and To = 170 MeV. 

These results were obtained in a stationary one- 
dimensional model with a single flow velocity. In real¬ 
ity different space-time sections of the overall FO layer 
are moving with respect to each other with considerable 
velocities, i.e. v « 0.2 — 0.7. Therefore, the superpo¬ 
sition of these parts of the FO layer wash out the very 
sharp peaks at small momenta, while the curvature at 
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higher momenta, although it is smaller, may persist even 
after superposition. There are several effects mentioned 
in the literature, which can cause such a curvature. The 
effects discussed in this section, arising from kinetic de¬ 
scription, may contribute to the curvature of the spectra, 
but we need a more realistic full scale, nonstationary 3- 
dimensional model to estimate the expected shape of the 
Pt spectra in measurements. Consequently, both the con¬ 
tributions of space-like and time-like sections of the FO 
layer have to contribute. 


0.1 X 1 X 10 X 



FIG. 15: The post FO distribution, in RFF. The 

calculation was made using the covariant escape rate, Wl for 
a finite, (L = 10A), FO layer. The subplots correspond to 
x = 0.1A, 1A, 10A respectively. The initial parameters are 
vo = 0 and To = 170 MeV. Contour lines are given at values 
represented on the figure. 


D. FFeeze out in layers of different thickness 

In this section we show the results of calculations per¬ 
formed with Wl escape rate, for different finite FO layer 
thicknesses. Some results of such an analysis have also 
been presented in 26]. 

In Figs. unnm we present the evolution of the tem¬ 
perature and flow velocity of the interacting component 
for L = 2A, 5A, 10A, 15A. We plot the resulting curves 
as function of x/L, what allows us to present them all 
in one figure. We clearly see, and this agrees also with 
our previous comparison to infinitely long FO, that by 
introducing and varying the thickness of the FO layer, 
we are strongly affecting the evolution of the interacting 
component. 

We can also study how fast the energy density of the 
interacting component is decreasing, see Fig. Since 
there is no expansion in our simple model, the evolution 
of the energy density is equivalent to the evolution of the 
total energy of the remaining interacting matter. We can 
see that the decrease of the energy density of the inter¬ 
acting component is exponentially fast, what justifies our 


way of getting approximate an solution for the interact¬ 
ing component, see section urn 



FIG. 16: The evolution of the temperature of the interacting 
component in RFF, for a baryonfree massless gas, calculated 
with the Wl escape rate for different FO layer thicknesses 
L = 2A, 5A, 10A, 15A. The initial parameters are t>o = 0.5 
and To = 170 MeV. 



FIG. 17: The evolution of the flow velocity of the interacting 
component in RFF, for a baryonfree massless gas, calculated 
with the Wl escape rate for different FO layer thicknesses 
L = 2A, 5A, 10A, 15A. The initial parameters are Vo = 0.5 
and To = 170 MeV. 


Figure HTH shows the final post FO transverse momen¬ 
tum distribution for different L. Despite the differences 
in the evolution of the interacting component, all the 
final post FO distributions look the same and are practi¬ 
cally indistinguishable. The difference between the result 
for a FO layer as thin as L = 2A and that for L —> oo 
limit shows up only in the low momentum region, and 
it is not significant enough to allow us to resolve layers 
of different thicknesses from experimental spectra. Thus, 
the thickness of the FO layer does not affect, as we have 
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FIG. 18: The evolution of the energy density of the inter¬ 
acting component in RFF, for a baryonfree massless gas, 
asB = 7 t 2 /10, calculated with the Wl escape rate for differ¬ 
ent FO layer thicknesses L = 2A, 5A, 10A, 15A. The initial 
parameters are no = 0.5 and To = 170 MeV. 


seen already in the previous section, the final post FO 
distribution, which is in fact the measured quantity! 


IV. CONCLUSIONS AND OUTLOOK 

In this paper we presented a simplified, but still non¬ 
trivial, Lorentz invariant freeze out model, which allows 
us to obtain analytical results in the case of a massless 
baryonfree gas. In addition the model realizes freeze out 
within a finite freeze out layer. 

We do not aim to apply directly the results presented 
here to the experimental heavy ion collision data, instead 
our purpose was to study qualitatively the basic features 
of the freeze out effect, and to demonstrate the applica¬ 
bility of this covariant formulation for FO in finite length. 

In Figs. 0 EH El and Figs. El ME! we compare re¬ 
sults, with the simple, cos 0 p , angular factor, and with 
the Lorentz invariant angular factor, ■ The differ¬ 

ences are insignificant, supporting our generalization. 

As it has been indicated in the previous publications 
E3, eh m, the final post FO distributions are non- 
equilibrated distributions, which deviate from thermal 
ones particularly in the low momentum region. The final 
spectra have a complicated form and were calculated here 
numerically. In large scale (e.g. 3-dim CFD) simulations 
for space-like FO the Cancelling Jiittner distribution 127^ . 
may be a satisfactory analytical approximation. 

Our analysis shows that by introducing and varying 
the thickness of the FO layer, we are strongly affecting 
the evolution of the interacting component, but the final 
post FO distributions, even for small thicknesses, e.g. 
L = 2A, look very close to our results for an infinitely 
long FO, first obtained in iiSiM. 

The results suggest that if the measured post FO spec¬ 
trum is curved, as shown in Fig. 1191 then it doesn’t 


matter how thick FO layer was, and we do not need to 
model the details of FO dynamics in simulations of colli¬ 
sions! Once we have a good parametrization of the post 
FO spectrum (asymmetric, non-thermal), it is enough to 
write down the conservation laws and non-decreasing en¬ 
tropy condition with this distribution function Q, (and 
probably with some volume scaling factor to effectively 
account for the expansion during FO). This Cooper-Frye 
type of description can be viewed from two sides. From 
experimental side, when we know the post FO spectra, 
we can extract information about the conditions in the 
interacting matter before FO. In theoretical, e.g. fluid 
dynamical, simulations such a procedure would allow us 
to calculate parameters of the final post FO distribu¬ 
tions to be compared with data. In this way our results 
may justify the use of FO hypersurface in hydrodynami- 
cal models for heavy ion collisions, but with proper non- 
thermal post FO distributions. 

At the same time, while the final distribution, f(p), is 
not sensitive to the kinetic evolution, other measurables, 
especially the two particle correlation function may be 
more sensitive to the details and extent of the FO pro¬ 
cess. 

The model can also be applied to FO across a layer 
with time-like normal. While several of the conclusions 
can be extended to the time-like case, it also requires still 
additional studies 0 - 

For realistic simulations of high energy heavy-ion reac¬ 
tions the full 3D description of expansion and FO of the 
system should be modeled simultaneously. We believe 
that our invariant escape rate, can be a basic ingredient 
of such models. 
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APPENDIX A 

The definition of K. n (z , w) function is: 

2” rd f°° i 

ICniz,w) = —z- n dxe~ x (x 2 - z 2 r~ * , (Al) 
( 2? U ! Jw 
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P x [GeV] 


FIG. 19: Final post FO transverse momentum (here p x ) distributions for the FO layers of different thicknesses. Calculations 
were done for baryonfree massless gas with escape rate Wl for L = 2A, 5A, 15A and with Woo (thick line). The initial conditions 
are the same as in Figs. nciirm To = 170 MeV, vo = 0.5. Distributions for the different FO layer thicknesses are very similar, 
with some difference in the low momenta region, which is shown in more detail in the ’’zoom in” subplot. The two thick lines 
correspond to some effective thermal distributions, with the corresponding parameters displayed in the plot legend. These are 
shown to illustrate the difference between obtained post FO distributions and thermal distributions. 


where in the case of w = z and n > — 1 the above formula 
will lead to the modified Bessel function of second kind, 
K n (z). Furthermore, the indefinite integral [28] is: 


I Z °' T(n,z)dz = 


z a T(n, z) — T(n + a, z) 


a 


, (A2) 


where r(ro, z) is the incomplete gamma function: 


F(n, z) = / dtt n 


(A3) 


The analytically not integrable functions G n (to) and 
H~{m) are defined as: 

G~{m) = j dpp [y/p 2 + m 2S j (A4) 


and 


i r 

( to ) = -— / 

n t > J'n+2 


dpp 


n+1 


(A5) 




Values of these functions for (n = 1,2) are given below: 

G,-(0) ^ ^ | <ipp 2 r(o,ip(l-i«)) (A6) 


= w^- iu) 


-3 


G 2 (°) = ^4 f Q d PP 3r (°. ■fPi 1 -J u j) ( A7 ) 


3 

2^ 


= 7TT ( X — j u ) ■ 


In the general calculation of the integrals in RFF, we 
change variables from p to z as given below: 


f‘OC _ 

/ dpf(p) e -l^P 2 +™ 2 -i 

Jo 

l d i' 

x / 7 T(juz- j\/z 2 - a 2 ) 


-jup) 


(AS) 


= JlT 


V z 2 — a 5 


JlT 


r*t 


u + j 


y/z 2 ~ a 2 ' 


x / 


7 T(juz + \J z 2 — a 2 ) 


In the massless limit, we have the G n (0) = H n (0). where z = 7 (y/p 2 + to 2 - jup)/T , a = m/T and b = 7 a. 
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